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I.  INTRODUCTION 


1.  STATEMENT  OF  THE  PROBLEM 

This  document  describes  the  theoretical  modeling  of 
turbulent  mixing  at  density  discontinuities  in  nonsteady 
compressible  flows.  It  is  well  known  that  if  a  density  dis¬ 
continuity  (or  a  strong  density  gradient)  occurs  in  a  pres¬ 
sure  gradient  of  the  opposite  sign,  then  the  flow  field  is 
hydrodynamically  unstable  in  the  Rayleigh-Taylor  sense 
(Ref.  1).  Small  perturbations  occurring  in  such  a  region 
will  amplify  with  time,  and  if  the  perturbations  become  large 
enough,  they  can  lead  to  a  local  breakdown  in  the  well-ordered 
flow;  i.e.,  they  can  lead  to  turbulence. 

One  of  the  more  interesting  examples  of  density  gradients 
working  against  pressure  gradients  to  cause  instabilities 

★ 

occurs  in  blast  waves  driven  by  solid  high  explosives  (HE). 
After  the  detonation  wave  breaks  out  of  the  charge,  the  HE 
combustion  gases  expand  to  a  high  velocity  (^6  km/s) ,  pushing 
an  air  shock  ahead.  One-dimensional  inviscid  calculations 
of  this  problem  (e.g..  Ref.  2)  show  that  a  positive  pressure 
gradient  is  formed  throughout  the  flow  field.  Such  calcula¬ 
tions  also  indicate  that  there  is  a  large  density  jump 
(p„_/p  ,  can  be  as  large  as  70)  across  the  contact  surface. 

HE  cl-liT 

However,  high-speed  photography  of  HE  experiments  shows  that 
this  contact  surface  (which  theoretically  should  be  smooth) 
actually  develops  an  irregular  shape  indicating  the  growth 
of  instabilities.  Evidence  of  mixing  at  the  contact  surface 
can  be  inferred  from  test  results  which  show  that  the  HE 
gases  react  with  the  shocked  ambient  gases  and  release  heat 


* 

Similar  gradients  can  occur  in  "shock-tube-type"  flows  driven  by 
gases  at  very  large  initial  densities  and  across  progressing 
flame  fronts  (i.e.,  deflagration  waves). 
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(i.e.,  afterburn)  if  the  ambient  gases  contain  oxygen 
(Refs.  3,4).  Since  this  heat  release  can  be  of  the  same 
order  as  the  energy  released  by  the  detonation  wave  (e.g., 

TNT  postdetonation  energy  release  is  about  2.5  times  the  deto¬ 
nation  energy) ,  it  will  affect  the  blast  wave  flow  field  and 
must  be  taken  into  account  for  accurate  numerical  simulations 
of  such  flows. 

This  report  describes  the  modeling  of  such  mixing 
phenomena  by  the  use  of  a  k-e  turbulence  model.  The  funda¬ 
mental  idea  is  to  generate  turbulent  kinetic  energy  according 
to  the  above-described  Rayleigh-Taylor  mechanism,  which  will 
then  produce  turbulent  mixing  in  the  region  of  the  contact 
surface.  Note  that  this  turbulence  mechanism  is  driven  at 
least  initially  by  pressure  and  density  gradients  normal  to 
the  contact  surface,  and  not  by  shear.  Turbulence  models 
based  on  such  normal  gradients  have  not  been  published  previ¬ 
ously.  This  report  describes  the  development  of  the  k-e 
model  and  presents  closures  for  the  terms  arising  in  such 
variable  density  flows. 

2.  COMMENTS  ON  FLOW  STABILITY 

Flow  instabilities  generate  and  sustain  turbulence 
(Ref.  5) .  The  following  elementary  argument  is  used  to 
identify  the  origin  of  an  instability  in  the  present  flow 
situation,  which  can  lead  to  the  outbreak  of  turbulence.  It 
should  be  stressed,  however,  that  these  arguments  are  of  a 
qualitative  nature  and  do  not  lead  in  any  way  to  a  quantita¬ 
tive  model  for  the  prediction  of  the  time-evolution  of  the 
flow. 

First,  a  definition  of  instability  is  loosely  given  here 
as  that  condition  under  which  a  perturbation  in  the  state  of 
a  medium  of  fluids  in  equilibrium  (or  quasi-equilibrium) 
results  in  a  departure  from  the  state  of  equilibrium  larger 


. .  I  'T  H  a 
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than  the  driving  perturbation.  Consider  now  two  fluids  of 
densities  pa  and  pg  (pa  >  pB )  in  a  state  of  quasi-equilibrium 
and  subject  to  a  pressure  gradient  in  the  direction  normal  to 
the  interface.  The  theorem  of  Bjerknes  (Ref.  5)  states  that 
the  rate  of  generation  of  circulation  r  (a  measure  of  vorti- 
city)  is  given  by 


dr  Cl 

■atm-jr  p'd-  (1) 

where  the  integral  is  taken  over  a  closed  contour.  Consider 
an  integration  path  1-2-3-4  in  Figure  1  that  encloses  a  finite 
portion  of  the  interface  separating  the  two  fluids,  and  where 
points  2  and  3  and  points  1  and  4  are  arbitrarily  close  but 
straddle  the  surface.  Eq.  (1)  then  yields 

=  ■  ^(P2“Pl}  -  |(P3"P2)  -  ^(P4“P3}  "  F(P1_P4)  (2) 

where  p  is  some  average  density  across  the  interface  and 
PA  >  p  *  PB*  S;*-nce  t*le  fluids  are  in  a  state  of  equilibrium 
and  there  is  no  pressure  gradient  tangential  to  the  interface, 
then 


Pi  =  P' 


=  P3  " 


and  Eq.  (2)  reads 


dr 

dt 


0 


Thus,  no  vorticity  is  generated  and  the  flow  retains  its 
equilibrium. 


Unperturbed  interface 


Figure  1.  Description  of  integration 

paths  of  Eq.  (1)  for  perturbed 
and  unperturbed  interfaces. 


Consider  now  a  perturbation  in  the  position  of  the  inter¬ 
face  such  that  points  2  and  3  are  removed  to  new  locations 
2*  and  3'.  A  pressure  gradient  tangential  to  the  interface 
is  now  set  up;  Eq.  (2)  then  reads 

dr  _  i  ,  .1,  .  i,  .1.  . 

dt  ‘  pA(p2,_Pl)  "  p(p3,_P2,)  "  pB  P4~P3 '  *  p(pl"p4) 


But  p1  =  p4  =  p2  and  p2 ,  =  p3<;  hence. 


For  the  case  when  p_  >  p_, 

A  B 


dr 

at  ~ 


const. 
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In  this  case  the  constant  is  positive,  and  3  p/'-  s  is  the 

pressure  gradient  normal  to  the  unperturbed  interface.  When 

3 p/3s  is  positive,  a  positive  circulation  is  generated  which 

tends  to  distort  the  interface  further — an  unstable  situation. 

Conversely,  when  3  p/9  s  is  negative,  a  negative  circulation  is 

generated  which  tends  to  restore  the  interface  towards  its 

equilibrium  state — a  stable  situation.  The  flow  being 

considered  in  the  present  report  corresponds  to  the  positive 

pressure  gradient  case  with  p,  >  p_,  and  hence  it  appears  to 

A  o 

contain  inherent  instabilities.  It  will  be  shown  later 
with  a  totally  different  approach  that  this  mechanism 
(a  pressure  gradient  acting  in  a  variable  density  flow  field) 
is  one  of  the  sources  of  turbulence  generation. 

3.  ON  THE  APPLICABILITY  OF  TURBULENCE  MODELING 

In  the  preceding  section  it  was  demonstrated  that  flow 
instabilities  can  arise  in  the  problem  under  consideration. 
What  is  not  clear,  however,  is  how  such  instabilities  evolve 
once  they  have  developed  and  whether  they  lead  to  the  tran¬ 
sition  to  full  turbulence  in  the  normally  accepted  sense 
(Ref.  6).  In  the  absence  of  detailed  experimental  data, 
doubts  must  arise  about  the  validity  of  classifying  the 
breakdown  phenomenon  at  the  interface  as  turbulence.  These 
doubts  are  enforced  when  it  is  remembered  that  the  macro¬ 
scopic  flow  field  evolves  so  rapidly  that  the  time  scale  of 
turbulence  fluctuations  becomes  similar  to  the  macroscopic 
time  scale.  The  applicability  of  turbulence  modeling  thus 
becomes  harder  to  justify. 

On  the  other  hand,  turbulence  can  be  viewed  as  a  random 
motion  subject  to  statistical  laws;  i.e.,  that  the  random 
motion  leads  to  rapid  mixing  and  property  exchange.  Since 
it  was  found  experimentally  that  very  rapid  mixing  occurs  in 
the  present  flow  case,  it  can  be  presumed  that  such  a  random 
statistical  motion  is  present.  Turbulence  modeling  endeavors 


to  simulate  the  exchange  processes  occurring  in  this  random 
motion  on  a  macroscopic  scale;  its  applicability  to  the  present 
case  can  therefore  be  justified  on  these  grounds.  However, 
this  argument  cannot  be  regarded  as  conclusive  until  it  is 
supported  either  by  experiment  or  by  analysis  of  the  micro¬ 
scopic  scale  motion — tasks  that  are  beyond  the  present  effort. 

4.  CHOICE  OF  TURBULENCE  MODEL 

It  is  assumed  here  that  the  reader  is  familiar  with 
turbulence  modeling;  those  who  are  not  can  find  exhaustive 
discussions  on  turbulence  and  its  modeling  in  may  texts,  such 
as  References  6  through  9. 

Most  turbulence  models  were  originally  developed  for  the 
special  case  of  constant  density  shear  flows.  A  comprehen¬ 
sive  review  of  these  models  and  their  closures  is  given  in 
Reference  9.  These  standard  models  contain  a  manifold  of 
fluctuation  correlation  terms  that  require  closure. 

Models  developed  for  compressible  shear  flows  are  based 
on  two  types  of  ensemble  averaging.  The  first  is  the  usual 
Reynolds  averaging,  which  results  in  many  terms  involving 
density  fluctuation,  each  requiring  closure.  Models  based 
on  this  type  of  averaging  will  inevitably  be  complex.  An 
example  of  a  model  that  utilizes  the  Reynolds  averaging  pro¬ 
cedure  can  be  found  in  Reference  10.  In  that  treatment,  many 
simplifications  and  restrictions  pertaining  to  the  atmospheric 
boundary  layer  were  introduced  at  the  outset  which  greatly 
simplified  the  model.  The  second  approach  is  the  mass- 
weighted  averaging  due  to  Favre  (Ref.  11).  It  results  in 
fewer  terms  at  the  expense  of  clarity  of  physical  interpre¬ 
tation  of  the  averaged  variables.  Examples  of  models  that 
are  based  on  mass-weighted  averaging  can  be  found  in  Refer¬ 
ences  12  and  13. 
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Although  the  k-e  model  (Refs.  7,14)  was  developed  for 
incorapressible  flow,  there  is  no  objection  in  principle  to 
its  generalization  to  flows  with  variable  density  as  was 
done  in  Reference  15  with  the  aid  of  mass-weighted  averaging. 
All  that  is  required  is  to  devise  model  closures  for  the  com¬ 
pressible  flow  terms.  The  model  in  its  incompressible  form 
has  been  extensively  used  over  a  very  wide  range  of  flow 
problems  and  has  proved  to  yield  reliable  predictions  in  many 
different  flow  situations.  It  was  applied  in  its  compres¬ 
sible  form  as  used  in  Reference  16  to  the  calculation  of 
supersonic  flow  with  imbedded  shocks.  With  all  these 

considerations  taken  into  account,  the  k-e  model  appears  the 
most  attractive  model  of  the  alternatives.* 

It  is  felt  that  none  of  the  existing  compressible  flow 
models  deals  adequately  with  the  "enthalpic"  generation  term 
of  turbulent  kinetic  energy,  which  is  regarded  in  the  present 
case  as  the  main  generator  of  turbulence.  In  these  models 
attention  has  been  focused  on  the  "kinetic"  generation  due  to 
stresses  and  strains  rather  than  the  enthalpic  generation  due 
to  density  variations  in  a  pressure  gradient  field  which  is  of 
major  relevance  here.  Therefore,  the  main  objective  of  this 
document  is  to  devise  a  suitable  closure  model  for  this  all- 
important  term. 

* 

More  sophisticated  models  which  solve  for  the  turbulent 
fluxes  were  not  seriously  considered  here  because  of  economic 
considerations.  Experience  shows  that  they  yield  moderate 
improvement  in  predictions  at  a  great  cost.  However,  they 
should  not  be  ruled  out  for  future  considerations. 
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II.  ANALYSIS 


1.  THERMODYNAMIC  RELATIONS 

The  following  analysis  includes  the  presence  of  up  to 
three  species  of  fluids;  the  generalization  to  a  system  con¬ 
sisting  of  more  species  is  a  straightforward  extension  of 
what  is  presented.  Different  fluids  will  be  labeled  with 
subscripts  A,  B,  and  C,  having  mass  mA*  mB,  and  mc,  respectively. 
The  mass  of  the  mixture  m  is  then  given  by 


m  =  mA  +  mB  +  mc 


(3) 


The  mixture  fractions  f  ,  f  and  f  are  defined  as 

A  o  C 


m 


f.  =  ~ 


f  =  — 
B 


A 

m 

mB 

m 


m 


=  — 


C 

m 


(4) 


From  Eqs.  (3)  and  (4)  it  is  clear  that 


f.  +  fn  + 
A  B 


fC  = 


(5) 


If  pA>  pB ,  and  pc  are  the  densities  of  the  species  A,  B,  and  C 
at  the  mixture  pressure  p  and  temperature  T  and  if  p  is  the 
density  of  the  mixture,  then  by  the  law  of  additive  volumes 
(Ref.  17),  one  obtains  the  following: 

V  =  VA(p,T)  +  VB(p,T)  +  Vc(p,T) 
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A 


iia— aaaii 


Hence, 


The  last  relation  can  be  rewritten  with  the  aid  of  Eq.  (4)  as 


Compressibility  Relations 

When  the  fluids  are  compressible,  the  following  relations 
may  be  invoked: 


P 


A 


P  ss  ...J?— 

B  RbT 
p  =  _E_ 

WC  RCT 


(7) 


where  it  has  been  assumed  that  the  fluids  act  as  perfect 
gases  with  gas  contants  R^,  Rg  ,  and  Rc-  The  assumption  of 
perfect  gas  behavior  is  not  a  prerequisite  and  may  be  replaced 
by  other  relations  if  required.  The  mixture  gas  contant  R  is 
given  by  the  mass-weighted  average 

R  '  £ARA  *  £BRB  *  £CRC  (8) 


Similar  relations  hold  for  the  specific  heat  quantities  C 
and  C^.  Combination  of  Eqs.  (6),  (7),  and  (8)  yields 


P 


(9) 
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2.  THE  TRANSPORT  EQUATIONS  IN  MASS- AVERAGED  FORM 


As  stated  earlier,  the  procedure  for  mass- averaging  was 
introduced  in  Reference  11  to  simplify  considerably  the 
resulting  equations  and  to  eliminate  a  great  many  of  the 
terms  involving  density-fluctuation  calculations.  The  basic 
definitions  are  now  presented,  and  the  transport  equations 
of  interest  are  then  stated. 

Let  the  instantaneous  value  of  a  dependent  variable  <j>  be 
composed  of  two  components:  a  mean  value  ^  (where  v  denotes 
mass-averaging)  and  a  fluctuating  component,  <j>'.  Thus, 


=  $  +  <*>' 


The  mass-averaged  quantity  $  is  defined  by 


(10) 


where  the  overbar  denotes  ensemble  averaging  (Refs.  11,15). 
<j>  may  stand  for  any  variable  except  p  and  p  for  which  the 
Reynolds  averages  are  used.  Thus, 


P  =  P  +  P’ 


and 


p  =  p  +  p 1 


It  follows  that 


P<)>'  =  0 


(11) 
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and  that 


P'Q' 

P 


t  0 


(12) 


The  transport  equations  for  mass,  momentum,  and  energy  of  the 
mixture  in  conjunction  with  the  transport  equations  for  the 
mass  fraction  can  be  written  in  mass-averaged  form  in 
Cartesian  tensor  notation  as  follows: 

Continuity 


+ 


0 


(13) 


Momentum 


■k(5“i>  +  fe-'evi’  -  srf?5ij  -  Tij  -  (14) 


where 


/9u.  3u.  ,  3uk\ 

Tij  5  +  3x.  -  3  6ij  3^/ 


and  u  is  the  molecular  viscosity. 
Energy 


-  a  (pfe+p,  pu'h') 

J  J  '  h  3  ' 


(15) 


'Xj 


'Xj 


ui  ("ij“  puiuj  *  +  ^f^j^PU^j) 


+  2f  (16) 
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where  l5  is  the  mass-averaged  total  energy,  h  is  the  specific 
enthalpy,  ah  is  the  Prandtl  number  for  molecular  heat  diffu¬ 
sion,  and  Q  is  the  external  source  of  energy  (heat  of 
reaction  when  present) .  The  instantaneous  total  energy  is 
given  by 


E  =  e  +  2  UiUi 


where  e  is  the  specific  internal  energy.  Taking  the  mass- 
weighted  ensemble  average  of  the  above  relation  gives 


*  -v.  1  ^  %  1  Puiui 

E  =  e  +  t  u.u.  +  ±  — ±-=- 

2  i  i  2  p 


where  the  last  term  is  the  average  kinetic  energy  of  turbulence, 
k.  The  last  expression  can  then  be  written  as 


1  V  7  Tj  <v, 

E  =  e  +  -=•  u .  u .  +k 


Mass  Fraction 


3  -c 

at  Pfa 


3  ^ 

3x~  puj  f a 


,3  \  u  3fa 

3Xj  |_aa  3x j 


pu  -  f ' 
D  a 


where  subscript  a  denotes  the  species  (A,  B  or  C) ,  a  is  the 

a 

Prandtl  number  for  the  molecular  diffusion  of  fluid  component 
a,  and  is  the  source  or  sink  of  component  a  (due  to  chemi¬ 
cal  reaction,  for  example). 

Also  needed  is  the  equation  for  the  transport  of 
turbulent  kinetic  energy  k,  which  is  as  follows: 


Turbulent  Kinetic  Ener 


_pk  +  3— (Pu.k)  =  Tx- \Ti jui  “  Puik  -  *1*1) 


_  a  a,  a  3u!  3u! 

PUlu'  5— -  U,  -  u!  -t—  p  -  p'^-i  -  T:.T-i 


i  j  3xj 


1  9xi 


ijSXj 


It  should  be  noted  that  the  preceding  transport  equations  are 
exact  insofar  as  the  original  instantaneous  equation  from 
which  they  are  derived  is  exact.  For  flows  without  fluctua¬ 
tions,  all  the  (')  quantities  vanish  and  the  equations  reduce 
to  their  laminar  forms.  The  significance  of  the  terms  involv¬ 
ing  fluctuation  correlations  in  Eqs.  (14)  to  (19)  is  discussed 
and  explained  in  the  many  existing  texts  on  turbulence  such  as 
References  6-9,  11,  and  15.  The  most  notable  of  the  correla¬ 
tion  terms  are  those  of  the  form  where  4>  may  stand  for 

any  of  the  dependent  variables.  These  terms  are  universally 
regarded  as  representing  the  turbulent  fluxes  of  property  t}>. 

The  aim  in  turbulence  modeling  is  to  provide  the  most  accurate 
and  universal  representations  of  these  terms,  which  is  called 
turbulence  model  closure. 

A  significant  difference  between  the  compressible  and  the 
incompressible  equations  is  the  presence  of  two  extra  terms, 
-u^(3p/3x^)  and  -p ' (3u^/9x^) ;  in  the  former  case,  both  terms  can 
be  regarded  as  turbulent  energy  generation  terms.  It  will 
be  seen  later  that  the  first  of  these  terms  plays  the  all- 
important  role  in  the  present  flow  problem. 

3.  THE  k-e  MODEL 

The  k-e  model  (Ref.  7)  employs  the  hypothesis  of  eddy 
diffusivity  for  the  evaluation  of  the  turbulent  exchange 
processes.  Thus,  for  any  turbulent  flux  -pu^p  for  scalar 
property  $,  the  expression  is  replaced  by 


(20) 


-Pu|<J> ' 


where  p  is  the  eddy  viscosity  to  be  defined  later  and  o 

t  t  t  (p 

is  the  turbulence  Prandtl  number  for  the  diffusion  of  scalar 
property  <(>.  For  the  turbulent  exchange  of  momentum,  the 
turbulent  flux  is  represented  by 


-puiuj  E  Vij  ‘  I  6ijkp  (21) 

where 

3u.  3u .  .  3uk 

£ij  “  3x .  +  3x.  3  (22) 

j  i  * 

Thus,  many  of  the  correlation  terms  are  replaced  by  terms 
involving  gradients  of  the  mean  quantities.  There  remain 
still  a  few  correlation  terras  that  require  closure;  they  are 
dealt  with  now,  one  by  one. 

The  term  -1/2 [uTTpuTUT) ]  in  the  energy  Eq.  (16)  is  almost 
the  same  as  -pu^k  which  is  as  already  discussed,  the  turbulent 
flux  of  k.  Hence,  the  term  may  be  replaced  by  the  following: 


(23) 


The  terms  3/3x^  “  u j_P '  )  Ec3-  (19)  are  normally  regarded 

as  diffusion  ones  and  are  lumped  with  -puTk.  Also  in  Eq.  (19) 
the  terms  -u|(3p/3x^)  -  p' (3u|/3x^)  require  closure.  These 
are  regarded  as  generation  terms  arising  from  the  variations 
in  density;  they  warrant  special  attention  for  this  problem 
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► 


and  are  therefore  dealt  with  in  a  separate  section.  For  the 
moment,  however,  they  will  be  lumped  with  the  kinetic  genera¬ 
tion  term  to  give  the  total  generation  term  G  as 


G 


3u. 

1 


3u ; 


-°uiuj  IxT  ~  ui  fx“  "  p'  3icT 


(24) 


Finally,  the  term  tTTTSulySxT ) i s  regarded  as  the  dissipation 
rate  per  unit  mass  of  the  turbulent  energy  and  is  represented 
by  a  new  variable  e  which  is  the  volumetric  dissipation  rate 
of  turbulent  energy.  Hence, 


3u! 

T  - i. 

ID  3Xj 


PE 


(25) 


The  variable  e  must  now  be  determined.  To  this  end,  a  trans¬ 
port  equation  for  e  is  formulated  and  solved  in  conjunction 
with  the  rest  of  the  transport  equations.  This  equation  takes 
a  form  similar  to  that  for  k  in  Eq.  (19) ,  and  is  given  by  Refer¬ 
ence  8. 


3_ 

3t 


pe 


(Pu^e) 


hr  far-  it)  *  ciGeA  +  c2DE2/k  (26) 


where  and  C2  are  empirical  constants. 

It  remains  to  determine  the  eddy  viscosity  ufc,  which  is 
evaluated  from 


-  „  2  . 
yt  =  cupk  z 


where 


is  another  empirical  constant. 


(27) 


The  model  is  now  closed  except  for  the  enthalpic 
generation  terms  in  Eq.  (24)  which  are  dealt  with  in  the 
following  section. 

4.  CLOSURE  FOR  ENTHALPIC  GENERATION 

The  generation  G  of  turbulent  kinetic  energy  is  given  by 
Eq.  (24)  as: 


G  = 


.a. 

3u.  -- 

-°ui"j  -sq  -  ui  sir 


The  first  term  in  the  expression  is  the  kinetic  production 
term  due  to  work  by  turbulent  stresses  acting  in  a  strain 
field;  the  closure  for  the  turbulent  stresses  is  given  by 
Eqs.  (21)  and  (22).  The  remaining  two  terms  in  the  expres¬ 
sion  for  G  are  the  enthalpic  generation  due  to  the  work  of 
fluctuating  mass-velocity  against  a  pressure  gradient  and  the 
work  of  fluctuating  pressure  against  a  fluctuating  strain 
field,  respectively. 

For  the  problem  under  consideration,  the  pressure  gradient 

is  considerable  and  the  density  fluctuations  are  very  large 

at  the  interface  because  of  the  high-density  gradient  across 

it.  The  term  -u| (3p/3x^)  is,  therefore,  expected  to  be  a 

major  contributor  to  the  generation  of  turbulence  at  the 

interface.  On  the  other  hand,  the  pressure  fluctuations  are 

not  thought  to  be  very  large,  and  smaller  still  will  be  the 

strain  fluctuations  3u!/3x..  The  term  p’ (3u!/3x. )  may  in 

ir  ii  * 

this  case  be  neglected  in  comparison  to  the  previous  term. 


*The  term  p' ( 3u! /3x . )  was  analyzed  as  a  diffusion  term  in 
Reference  18  an&  noi  as  a  generation  term.  It  is  not  clear 
how  this  conclusion  was  reached  as  no  justification  was  given. 
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The  task  now  is  to  provide  a  closure  for  the  term  -uT(3p/3x^), 
which  may  be  rewritten  with  the  aid  of  relation  Eq.  (12)  as 
(pTuT/p) (3p/3x^).  In  what  follows,  two  independent  paths  are 
taken  to  obtain  such  closure;  it  will  be  shown  that  the  result¬ 
ing  closures  are  almost  equivalent. 

a.  First  approach — To  arrive  at  an  expression  for  the 
correlation  p"ruT,  one  may  use  the  arguments  presented  in 
Reference  19.  Since  mass-weighted  averaging  is  employed 
presently,  the  derivation  below  is  a  slight  modification  of 
the  one  presented  in  that  reference. 

The  density  fluctuation  may  be  written  in  terms  of  mass- 
weighted  quantities  as 


P' 


L_  3P 

Sk  p  9x 


(28) 


where  L  is  a  length  scale.  When  Eq.  (28)  is  multiplied  by 
u!  and  then  ensemble-averaged,  the  following  relation  is 
obtained: 


P  1  u!  = 
x 


r  pu!u'. 

L  i  j  3p 


/k 


3Xi 


(29) 


Now  L  is  given  in  Reference  20  as 


L  =  C 


3  2 


/k 


(30) 


and  -puTuT  is  here  modeled  by  Eq.  (21).  Hence, 


9p 


*  e.. 


p '  u!  =  -C,  5.  .  -p-  +  ^  C,  e.  .  - — 

l  3  p  13  dx.  2  3  -2  13  3x. 

3  p  k  3 
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where  is  another  empirical  constant.  The  enthalpic  gener¬ 
ation  term  thus  becomes 


_^r  3fL  =  ^  f.  ]££ 

1  8xi  p2  L  ij  2  ijJ  3xj  Sxi 


When  3p/3x^  is  negative  and  3p/3x.  is  positive,  as  is  the 
present  case,  for  i=j  the  term  in  the  square  bracket  is  almost 
certain  to  be  positive  since  the  normal  strain  field  is  small; 
hence,  the  generation  term  gives  a  positive  contribution.  For 
the  case  when  the  sign  of  either  3j5/3x.  or  3p/3x.  is  reversed, 
the  term  reverses  its  role  and  acts  to  destroy  turbulence 
energy.  The  action  of  the  term  does,  therefore,  follow  the 
trends  already  discussed  in  Section  1-2. 

b.  Second  approach — In  the  second  approach,  the  density 
fluctuation  is  obtained  from  an  equation  that  relates  the 
density  to  other  dependent  variables;  Eq.  (6)  serves  this 
purpose,  a  displacement  p'  is  considered  to  result  in  a  dis¬ 
placement  in  all  other  variables  such  that  the  ensemble 
average  of  fluctuations  in  mass-weighted  terms  vanishes.  Thus 
p'  results  in  displacement  pf^/p  and  PT'/P  and  so  on.  Eq.  (6) 
then  gives 


It  should  be  noted  that  the  first  term  on  the  right-hand 
side  of  Eq.  (32)  is  due  to  the  variations  in  density  regard¬ 
less  of  the  effects  of  compressibility;  it  is  present  even 
if  the  flow  is  incompressible.  The  second  term  expresses 
the  effects  of  compressibility  of  the  different  gas  components 
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From  Eq.  (7)  the  fluctuations  p',  p'  and  p'  may  be  obtained 

Ad  L. 

in  terms  of  fluctuations  in  T  and  p.  Thus, 


Pg  _  p'  pT ' 
Pa  P  PT 


(33) 


where  a  stands  for  A,  B  and  C. 

Substitution  of  Eq.  (33)  into  Eq.  (32)  gives 


Combination  of  the  last  equation  with  Eq.  (6)  then  gives 


(34) 


Eq.  (34)  is  now  used  to  replace  p'  in  the  enthalpic  production 
term  which  is  restated  here  as 


Enthalpic  Generation 


_sr  |e_  , 

i  9x . 
i 


Introduction  of  Eq.  (34)  yields 


-st|e_  = 

i  9x. 


i 


r 
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Thus  far,  the  treatment  of  the  effects  of  density 
fluctuations  has  been  exact.  But  now,  approximations  are 
made  (or  modelled)  in  the  individual  terms  in  the  last 
expression.  According  to  the  eddy  diffusivity  model  that 
is  being  used  presently,  Eq.  (20)  can  be  used  to  replace 
many  of  the  terms  in  the  expression.  Hence, 


3P 

3x. 

1 


3  ft. 


B 


-  3x. 

PB  1 


p'  u! 

r  l 


P 


+ 


a  ,  —y 
t ,h  ph 


where  in  the  last  expression,  the  temperature  has  been 
replaced  by  the  specific  enthalpy.  Earlier,  it  was  argued 
that  pressure  fluctuations  are  negligible  compared  to  the 
rest  of  the  fluctuations;  therefore,  the  term  p' u!  can  be 
dropped.  The  final  expression  then  becomes 


iE_  = 

ui  3xi 


Mt  /  1  ^A  +  J_^B  + 

“t.fk  0B 


3ft 

3x . 


x 


(35) 


So  far,  two  different  expressions,  Eqs.  (31)  and  (35), 
have  been  derived  to  represent  the  same  term.  It  will  be 
shown  now  that  they  are  essentially  equivalent.  Consider 
Eq.  (31);  if  the  term  pertaining  to  the  rate  of  strain 
neglected,  the  expression  becomes 


is 


C3pt 

-2 

P 


3  P 
9xi 


3x. 

l 
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Differentiation  and  manipulation  of  Eqs.  (6)  and  (7)  give 


1 

-2 

P 


Therefore, 


_L 

P*  T*! 


■“i  Hr  -  Vt 


_1_  _1_  J_  3£_  _1 _ 

-  3x..  -  3x_.  -  3x_.  ~  --  9x4  3x. 


3h 


B 


PP 


ph 


The  last  expression  is  almost  identical  to  Eq.  (35);  the 
extra  term  involving  3p/3x^  is  not  obtained  in  Eq.  (35) 
because  of  the  neglect  of  the  term  p' u!  in  its  derivation. 

It  is  clear  that  the  two  different  approaches  lead  to  almost 
identical  expressions,  which  enhances  the  confidence  in  their 
correctness. 

A  final  recommendation  as  to  what  form  of  the  enthalpic 
generation  term  should  assume  may  be  put  forward.  Since  the 
neglect  of  term  e„  in  Eq.  (31)  leads  to  an  expression  very 
similar  to  that  in  Eq.  (35),  hence,  it  may  be  neglected, 
leaving  the  expression 


1E_  =  - 
3x.  -2 

i  P 


c3pt  3p_ 


3x .  3x . 


(36) 


As  Eq.  (36)  is  simpler  than  Eq.  (35),  the  former  is  recommended 
for  use  in  the  k-equation. 

5.  FINAL  SET  OF  EQUATIONS 

The  final  set  of  equations  will  now  be  summarized.  Since 
molecular  diffusion  plays  a  trivial  role  in  the  large-scale 
mixing  of  the  present  problem,  the  terms  expressing  its  effects 


a>|  qj 


have  been  dropped.  The  overbars  and  ^  will  be  suppressed, 
because  all  the  dependent  variables  are  taken  to  be  the  mean 
value  of  the  property  solved  for. 

Continuity 


ft  +  fxT  <»“!>  -  0 


(37) 


Momentum 


ft  pv  hq  foujui’ 


%  +  fj-  "t£ij  f  fx~  pk  ,38) 


where 


3u .  3u.  .  3u 

—  =  i  i  ^  ^  c  n 

ij  ~  3x.  ~  I  lj  3xn 


(39) 


Energy 


|_0E.|_OUj(E+P/p) 


3x  . 
3 


+  Q 
(40) 


where 


EBe  +  tt-u.u.  +  k 

2  ii 


(41) 
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Mass  Fraction 


3  *  .  3  ,  *  ,  _  3  /  ut  9fct\  .  „ 

3t  pfa  pu j f a  ~  3x^  | afc  f  3x .  f  Sa 


Turbulent  Kinetic  Enerc 


3  v  ^  3  /  3  /  Mt  3k  \ 

3t  P  3*j  j  3*j  U.k  3*:i) 


+  G  -  pe 


where 


3u.  C_ y 


{  2  ,  .  \  i  3  t  3p  3g_ 

K£ij  -  T  6ijk0)  9x7  7T"  '9x7  9x7 


Dissipation  of  Turbulent  Kinetic  Enerc 


|_  0S  +  |_  (0u.£) 


Turbulent  Viscosity 


_  /  yt  1§_  \ 

cj  K,s  9xj  1 


+  C1G  ^  "  C2P  V  U5) 


_ 

"t  ‘  V  T 


This  closes  the  set  of  equations. 

Turbulence  model  constants  C^,  C2,  ,  and  the  turbulent 

Prandtl  numbers  afc  etc. ,  have  been  derived  for  turbu¬ 

lent  shear  flows.  Values  for  these  constants  that  are  normally 
recommended  for  such  shear  flows  are  given  in  Table  1.  These 
values  can  be  used  as  initial  guesses  for  the  present  problem . 
However,  one  should  be  prepared  to  experiment  with  these 
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values  and  to  verify  their  validity  for  such  Rayleigh-Taylor 
mixing  problems  by  comparison  with  experimental  data.  A  new 
constant  has  been  introduced  in  the  present  document?  it 
is  roughly  inversely  proportional  to  the  turbulence  Prandtl 
number  (compare  Eqs.  (36)  and  (37)).  It  should  therefore  lie 
between  0.7  to  1.5. 


TABLE  1.  VALUES  OF  CONSTANTS,  AS  DERIVED 
FOR  SHEAR  DOMINATED  FLOWS 


Constant 

Value 

C1 

1.44 

C2 

1.92 

C3 

0.7-1. 5 

c 

.09 

u 

°t,h 

1.0 

G.  . 

1.0 

t,k 

Q 

rt 

Hi 

1.0 

a. 

1.3 

t,  £ 
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APPENDIX  A 

EXPANSION  OF  THE  GOVERNING  EQUATIONS 
FOR  SPHERICALLY  SYMMETRICAL  FLOW 


Eqs.  (37)  to  (44)  are  expanded  in  this  appendix  in 
spherical  coordinates  assuming  spherical  symmetry  of  the  flow. 

Continuity 
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Mass  Fraction 
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Turbulent  Kinetic  Energy 


(pit,  ♦  jj  (r2puk)  =  4  |r  (r2  0(;‘k  ,*) 
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Dissipation  of  Turbulent  Kinetic  Energy 
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The  constants  ,  C 2,  etc.,  were  given  in  Table  1. 
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